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Abstract. Reliable and exact assessment of visibility is essential for safe 
air traffic. In order to overcome the drawbacks of the currently subjective 
reports from human observers, we present an approach to automatically 
derive visibility measures by means of image processing. It first exploits 
image based estimation of the atmospheric transmission describing the 
portion of the light that is not scattered by atmospheric phenomena (e.g., 
haze, fog, smoke) and reaches the camera. Once the atmospheric trans¬ 
mission is estimated, a 3D representation of the vicinity (digital surface 
model: DMS) is used to compute depth measurements for the haze-free 
pixels and then derive a global visibility estimation for the airport. Results 
on foggy images demonstrate the validity of the proposed method. 


1 Motivation 

In order to guarantee for safe air traffic controllers have to rely on precise fore¬ 
casts and measurements of the current weather situation. Those have to be 
compiled and reported according to official regulations (e.g. visual flight rules 
(VFR) flight0) every 30 minutes in an international standard format [^] The 
decisions the controllers take do not only impact the security but also influence 
the air traffic and its economic repercussions in case of delays [T. 

All major airports operate dedicated sensor systems to assess the current 
weather situation. Besides “classic” parameters like wind, pressure, humid¬ 
ity, and temperature, there are point-like measurements of visibility and cloud 
cover information. For reporting the prevailing visibility in the airport vicinity, 
observers currently compile their reports based on integrating sensor measure¬ 
ments with visual observation of known landmarks, like buildings, mountain 
tops, etc. (Figure [l]). Human estimation is naturally very subjective to the 
individual observer and error-prone. Visibility sensors employ the principles of 
forward scattering to provide accurate and reliable information regarding vis¬ 
ibility changes in the atmosphere. These sensors (usually located at the start 
and end of the runways) are very precise but only give very local (e.g. at the 
position of the sensor) measurements. 

1 Visual meteorological conditions: http://en.wikipedia.org/wiki/Visual_ 

meteorological_conditions 

■^METAR: METeorologique Aviation Reguliere 
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Figure 1: Landmarks around the airport of Graz (Austria). 


Contrarily to Andreu et al. [2], the approach presented in this paper is not 
aiming at emulating the controller procedure which actually results in a limited 
number of visibility measurements (up to the number of defined landmarks 
scattered around the airport) but instead intends to get a dense (e.g. for each 
pixel in the image) visibility estimation. In order to quantify visibility, we 
propose a three steps procedure. First, we use a single image haze removal 
technique in order to recover the atmospheric transmission at every pixel of 
the image. Second, with the help of a very precise (e.g. lm resolution) digital 
surface model of the vicinity of the airport, a depth map is generated for the 
whole image. Finally, statistics about the depth of the scene points whose 
corresponding pixels are deemed visible gives the visibility in the field of view 
of the camera. 

A short overview of the relevant literature in this field is presented in Section 
[2| which is followed by the methods applied for estimating the atmospheric 
transmission in Section [3j Then, it is shown in Section [4] how, with the help of 
a very precise digital surface model of the airport vicinity, a depth map can be 
generated. In Section [5] we combine the estimated atmospheric transmission and 
the depth map to derive a global visibility estimation for the airport. Finally, 
Section [6] concludes with final remarks. 

2 Related work 

Weather and other atmospheric phenomena, such as haze, fog, mist, rain and 
snow, greatly reduce the visibility of distant regions in images of outdoor scenes. 
Quantifying the amount of atmospheric scattering or removing it (e.g. dehaz¬ 
ing or de-weathering), is a challenging problem, because the degree to which it 
effects each pixel depends on the depth of its corresponding scene point. This 
is leading to an under-constrained problem if the input is only a single image. 
Therefore, many methods propose using multiple images or additional informa¬ 
tion in order to alleviate that inherent problem. 

Recent dehazing methods like [3] and m are able to dehaze single images by 
making assumptions about the contrast difference in haze-free and hazy scenes 
as well as assumptions about the transmission and the surface shading. Even 
though those methods result in visually appealing images they are either not 
physically valid or fail in handling images with heavy haze. 

In order to remove haze, Schechner et al. Q2I successfully made use of 
multiple images taken with different polarizer orientations, which is difficult to 
adapt for standard webcam single images. 
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Depth-based methods, where depth information or heuristics about the scene 
are provided either by interactive user inputs or by photo realistic 3D models, 
have been successful in de-weathering images. Narasimhan et a 1. m addressed 
the question of de-weathering a single image using simple additional information 
(i.e. sky, vanishing point) provided interactively by the user. Kopf et a 1. [8] 
take advantage of the availability of accurate 3D models and try to estimate 
stable values for the haze function directly from the relationship between the 
colors in the image and those of the rendered 3D model. However, in our case, 
even if we do have access to a pretty precise digital surface model generated by 
laser scanning, this 3D model has no color information thus preventing the use 
of Kopf et a 1. ’s approach. 


3 Atmospheric Transmission 

3.1 Dark Channel Prior 

In computer vision, the model widely used to describe the formation of an image 
perturbed by haze is m, uni, mi: 

I (x) = J ( x ) t(x) + A (1 — t(x)) (1) 

where / is the observed intensity, J is the scene radiance, A is the global atmo¬ 
spheric light, and t is the atmospheric transmission describing the portion of the 
light that is not scattered and reaches the camera. The first term J(x)t{x) on 
the right-hand side is called direct attenuation, and the second term A( 1 — t(x)) 
is called air-light. When the atmosphere is homogenous, the transmission t 
is attenuated exponentially with the scene depth and can can be expressed as 
t{x) = where (3 is the scattering coefficient of the atmosphere and d the 

scene depth. 

He et a 1. [5 propose using a dark channel prior for single image haze re¬ 
moval. This prior is partially inspired by the well-known dark-object subtraction 
technique [3] used in multispectral remote sensing systems. The idea is that in 
most of the local regions which do not cover the sky, some pixels (called dark 
pixels) normally have very low intensity in at least one color (RGB) channel. 
In hazy images, the intensity of these dark pixels in that channel is mainly 
contributed by the air-light. Therefore, these dark pixels can directly provide 
an accurate estimation of the haze transmission. For an arbitrary image /, its 

dark channel I dark is given by I dark (#) = min ( min I c (y ) ) where I c is 

yeO(x) \ce{r,g,b} J 

a colour channel of / and is a local patch centred at x. Using the concept of 
a dark channel, if I is an outdoor haze-free image, except for the sky region, 
the intensity of F s dark channel is low and tends to be zero. Figure [2] shows a 
foggy image of the airport surrounding (left) and its corresponding dark channel 
computed over local patches of size 15x15 pixels (center). 


3.2 Estimating the transmission 


As in [5], A is automatically estimated from the most haze-opaque pixels (the 
top 0.1% brightest pixels in the dark channel). If the transmission t{x) can 
be assumed constant (denoted as t(x) ) over the local patch U(x), Eq. ([!]), 
normalized by A for each colour channel, could be rewritten as: 


min 

yett(x) 


( min 

\c£{r,g,b} 



t(x) min 

y£tt(x) 


( min 

\c£{r,g,b} 



+ 1 - t(x) 


( 2 ) 
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Figure 2: Foggy image, its dark channel computed over 15x15 patches, and 
transmission map after guided filtering with the source image. 


By definition the scene radiance J should be haze-free, its dark channel is close 
to zero so Eq. © leads to: 


t(x) = 1 — min 

y£tt(x) 


( min 

\c£{r,g,b} 


nv) 

A c 



( 3 ) 


where min ( min 1 & + 1 ) is the dark channel of the normalized hazy 

yefl(x) \ce{r,g,b} J 

image. In practice, even on clear days the atmosphere is not absolutely free 
of any particle meaning that haze still exists when we look at distant objects. 
So, a very small amount of haze could be optionally kept for distant objects by 
introducing a constant parameter cj(0 < uj < 1) into Eq. & 


t(x) = 1 — oj min 

y£tt (x) 


I min 

\c£{r,g,b} 


nv ) 

A c 



( 4 ) 


Note that this method for estimating the atmospheric transmission is based on a 
prior (Dark Channel) reflecting the statistics of outdoor images. That prior may 
become invalid on objects similar to the atmospheric light (i.e. white walls) or 
having no shadow cast on them. As a result, this method might underestimate 
the transmission of these objects and overestimate the haze layer. However, 
images of the vicinity of an airport do not contain many of such problematic 
objects and have plenty enough shadow casting. 


3.3 Refining the transmission 

Using Eq. © for computing the transmission results in a very coarse map. Since 
the transmission is not always constant over local Q(x) patches, the map has 
to be refined in order to capture the sharp edge discontinuities and outline the 
profile of the objects. We follow the approach of He et a 1. [5] applying a guided 
filter (with the original hazy image used as the guidance image) to refine the 
transmission map |B]. Figure [2] (right) shows the corresponding refined version 
of the transmission map. 

4 Depth map generation 

For this study we identified 16 ground control points (GCPs) manually for a part 
(i.e. camera field of view) of the airport vicinity (Figure [ 3 ]). The 3D reference 
coordinates of the GCPs were derived from a combined LiDAlOSRT]V0DSM. 


3 Light Detection and Ranging 
4 Shuttle Radar Topography Mission 
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In the northern part in the vicinity of the tower of the airport Graz we therefore 
used a LiDAR DSM which was acquired with 4pts per sqm in the year 2009. For 
more distant areas we used the globally available version 4 of the SRTM DSM, 
which was acquired using single pass interferometry at C-band with a ground 
resolution of about 90m in the year 2000 [7]. Figure [3] shows a painted relief of 
the generated combined LiDAR SRTM DSM. To facilitate the identification of 
the buildings given in the LiDAR DSM we additionally used Google map (See 
Figure [4]) and here with especially oblique views. 



Figure 3: Painted relief of the combined LiDAR SRTM DSM south to the tower 
(position outlined at the center of the red circle in the lower left part of the 
image) of the airport of Graz (Austria). Blue areas indicate no data. 



Figure 4: Identification of GCPs in one image of the webcam of the airport 
Graz (Austria) (left) and combined LiDAR-SRTM and Google map (right). 

Based on the well-known (linearized) co-linearity equations we used the GCP 
information to orient this image . Besides the extrinsic orientation parameters 
(camera position and rotation) we limited the intrinsic parameters to the focal 
length only. For a more precise distance map also the principal point and camera 
distortion parameters could be derived. A coarse position of the webcam and 
the knowledge that the camera is pointing to the South were sufficient to set-up 
the orientation with a root mean square error of less than 4 pixels (with respect 
to the GCPs reprojection error). Once the image is oriented, for every pixel 
of the image the first intersection of the line-of-sight (LOS) with the combined 
LiDAR-SRTM DSM is calculated and the distance assigned to the raw distance 
map. More precise, we restrict ourselves only to image pixels below the detected 
horizon. As the combined LiDAR-SRTM DSM is only 2.5 dimensional this 
calculation fails if the LOS is below e.g. the roof of a building. In this case the 
intersection is found with an object behind the actual building. We therefore 
applied a simple reconstruction of facades by projecting the eaves line to the 
ground. Then a simple search in the image line direction ensures that the range 
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direction is always decreasing (when starting from the horizon). Figure [5] shows 
the raw distance map of the airport Graz (Austria) and corrected distance map. 
All distances are given in logarithmic scale ranging from less than 60 m to more 
than 16 km. 



Figure 5: Raw depth map of the airport Graz (Austria) (left) and corrected 
depth map (right). Depths are given in logarithmic scale from 60m up to 16km. 


5 Visibility estimation 

Figure [6] shows the result of detecting haze on two hazy images taken at different 
day time under different weather conditions (low-level clouds and fog banks). 
The border of the haze is delineated by finding the boundaries of the region 
obtained by thresholding the transmission map computed using the methods 
described in Section [3] 

As can be seen on the both examples of Figure [6j the border between hazy 
and non hazy regions of the image is delineated by using a threshold of 0.75 
(heuristically estimated) on the corresponding transmission maps. On the right 
hand side, one can notice that the fog bank over the runway is accurately seg¬ 
mented which is of importance for air traffic. Derivating a visibility distance 



Figure 6: Highlighted (in red) haze border by using a global threshold of 0.75 
on the corresponding transmission maps. 

estimation out of the depth and transmission map is done by statistic analysis 
on the pixels of the depth map whose transmission value in the corresponding 
transmission map is above a given threshold value. One can then either choose 
as visibility distance the overall depth maximum of the concerned pixels or, in 
order to eventually avoid outliers, a simple percentile statistics over the gathered 
depths (i.e. the top 1%). 
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Table [T| reports visibility distance estimations from an initial evaluation to¬ 
gether with the Austro Control observers. It is done with a global threshold 
value of 0.75 on the transmission map and the distance corresponding to the 
99th percentile rank of the depth map histogram (with a bin distance of 10m) 
as resulting visibility distance. Globally, the results over a broad range of dis- 


Measured 

Reported 

(METAR) 

Measured 

Reported 

(METAR) 

Measured 

Reported 

(METAR) 

102m 

250m 

731m 

500m 

6768m 

6000m 

110m 

300m 

443m 

600m 

10537m 

10km 

371m 

300m 

791m 

800m 

12178m 

10km 


Table 1: Comparison of measured to officially reported visibility (METAR). 

tances (from as close as 250m up to 10km) are mimicking official reports well. 
Some discrepancies could be explained by the difficulty for the controller to 
evaluate the visibility distance due to the low number of landmarks in some 
distance ranges (e.g. in the range [600m, 800m], cf. Figure [l]) and also due to 
the subjectivity of the concept of visibility as human beings are defining it (e.g. 
when is landmark visible and when not?). Clearly, the value heuristically cho¬ 
sen as global threshold on the transmission has a major impact on the visibility 
results. No criticality analysis has been done on this parameter yet since the 
results presented here are mostly preliminary and mainly thought as a proof of 
concept. 

6 Conclusion 

In this paper, we have proposed an approach to automatically derive visibil¬ 
ity measures at airport sites by combining an image processing method with 
the digital surface model of the airport vicinity. The image processing method 
equivalent to an image haze removal technique enables us to first recover the 
atmospheric transmission at every pixel in the image. Then, with the help of a 
very precise digital surface model of the airport vicinity, a depth map is gener¬ 
ated for every pixel. Finally, gathering the depth of the pixels which are visible 
gives the visibility in the field of view of the camera. It was shown that this ap¬ 
proach gives preliminary results comparable to the visibility estimations given 
by human controllers. Further work has to be done in order to automatically get 
an optimal value of the transmission threshold (i.e. by exploiting local features 
for threshold derivation) which is a critical parameter of the proposed approach. 
Also, instead of being global, we will investigate threshold values that could be 
locally adapted to different parts of the image. 
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